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Abstract. Performance of the two commonly used numerical procedures, one based on artificial 
compressibility method and the other pressure projection method, are compared. These 
formulations are selected primarily because they are designed for three-dimensional applications. 
The computational procedures are compared by obtaining steady state solutions of a wake 
vortex and unsteady solutions of a curved duct flow. For steady computations, artificial 
compressibility was very efficient in terms of computing time and robustness. For an unsteady 
flow which requires small physical time step, pressure projection method was found to be 
computationally more efficient than an artificial compressibility method. This comparison is 
intended to give some basis for selecting a method or a flow solution code for large three- 
dimensional applications where computing resources become a critical issue. 


1 Introduction 

Various procedures can be selected for simulating incompressible flows depending on 
the choice of formulations, variables, discretization and iterative schemes. In the 
present study, performance of the two procedures commonly used for three - 
dimensional applications are compared: the first one based on an artificial 
compressibility method and the second one on a pressure projection method. 

The artificial compressibility method takes advantage of the advances made in 
conjunction with compressible flow computations. This approach relaxes the 
requirement of enforcing mass conservation equation rigorously at each time iteration, 
however, at the expense of introducing an artificial wave phenomenon. This approach 
can be viewed as a special case of a preconditioned compressible flow formulation. 
However, the computational efficiency is in general better than that of compressible 
flow solvers at the incompressible limit. This approach has been shown to be very 
robust in a wide range of applications [Kwak et. al (1986), Rogers et. al (1991)]. 

The first primitive variable method for incompressible flow was developed by Harlow 
and Welch (1965) using pressure projection. Numerous variants have been 
developed since. In this method, the pressure is used as a mapping parameter to 
satisfy the continuity equation. The usual computational procedure involves choosing 
the pressure field at the current time step such that continuity is satisfied at the next 
time step. The time step is advanced in multiple steps (fractional step) which is 
computationally convenient. However, governing equations are not coupled as in an 
artificial compressibility approach. This will affect the robustness and limit the 
maximum allowable time step size. Since this approach is time accurate, there are 
cases where the fractional step solver is computationally more efficient compared to 
the artificial compressibility method [Rosenfeld et. al (1991), Kiris and Kwak (1996)]. 

Various numerical algorithms associated with these methods have been developed 
along with accompanying flow solvers. In three-dimensional applications computing 
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time requirement is one of the key issues in selecting flow solvers. In the present study 
it is intended to give some bench mark comparisons on computing efficiencies 
between the two methods discussed above. 


2. Artificial Compressibility vs. Pressure Projection Formulations 


Governing equations for these two methods are presented here to facilitate the 
discussion of results. The artificial compressibility method by Chorin (1967) can be 
written as 


J _^ + ^L = Q 

(5 dt dx i 


du, 


dt 


i- = -Vp + h i 


( 1 ) 

(2) 


where usual conventions are used with representing convection and viscous terms. 
At steady state the pressure term in equation (1) drops out and thus incompressibility 
is recovered. For time accurate computations, this has to be repeated at each time 
level to maintain incompressibility at each time step. 


The pressure projection procedure on the other hand requires solving the following 
pressure Poisson equation to satisfy the continuity at the next time level: 

(3) 

Numerically intermediate velocity field is computed first then pressure correction is 
computed using the pressure Poisson equation. Once the pressure correction is 
computed, new pressure and velocities are calculated. In this approach the Poisson 
solver portion is usually the most expensive part of the computational procedure. 
Therefore, accelerating convergence for Poisson solver has been the focal point of 
many previous studies. 


In the present paper, the two methods are compared using solvers developed at 
NASA Ames Research Center. The artificial compressibility method is represented by 
INS3D-UP code which utilizes line relaxation in conjunction with upwind differencing. 
The pressure projection method is represented by INS3D-FS which utilizes Runge- 
Kutta time stepping in conjunction with finite volume scheme on staggered grid 
arrangement. Code performance is a function of many factors such as iterative 
scheme, spatial differencing, coding and choice of parameters like artificial 
compressibility. However, the performance figures presented next will provide some 
indications on what one might expect from similar codes in computing steady-state or 
time-accurate problems encountered in engineering practices. 


3. Computed Results 

Two test problems are selected: a wake vortex problem for steady-state solutions, and 
a curved duct flow with pulsatile inflow condition. The wake vortex problem is selected 
because it represents aerodynamic and hydrodynamic problems involving wingtip. 
The high gradient in the flow field poses significant challenges for computing this flow. 
The pulsatile flow in curved duct is chosen to represent internal flows frequently 
encountered in propulsion and bioengineering. 
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In Figure 1 , tip vortex generated from a wing is shown from an earlier numerical and 
experimental study by Dacles-Mariani et. al (1993). A round tip NACA 0012 wing is 
placed at 10 degrees angle of attack at /?, = 4.6xl0 6 . The present computational 
domain starts at the trailing edge of this wing. The two codes are compared for 
obtaining steady-state solution using 35x81x81 grid. The flow field has a highly 
concentrated vortical core region which poses challenges in implementing the right 
turbulence model and grid distribution. In Figure 2, convergence history using INS3D- 
UP and INS3D-FS is shown. For steady-state solutions, INS3D-UP code required 2 
hours of Cray C-90 CPU time while INS3D-FS needed 4.2 hours to lower the residual 
down to 10" !3 . The pressure projection required smaller time step to maintain the 
projection accuracy which resulted in approximately twice more iterations than the 
artificial compressibility method. In Figure 3, results using various production terms in 
the Baldwin-Barth turbulence model (denoted as BB) are compared using INS3D-FS. 
Also shown is a coarse-grid solution using lower order spatial differencing. The 
corresponding velocity profiles at x/c=1 .447 are shown in Figure 4. Figure 4 shows that 
the large differences among different experiments in Figure 3 are mainly in the core 
region, and in the outer region, the computed velocity compares well with experiment. 
This study shows the relative importance of selecting the right turbulence model as 
well as grid resolution in computing tip vortex flow . In Figure 5, the velocity 
magnitudes from the two codes and experiment are plotted at x/c=1 .587. Full details of 
this study can be found in Dacles-Mariani et. al (1993) and Kiris and Kwak (1996). 

For an unsteady test, a square duct with 120 degree bend was chosen. The inflow as 
shown in Figure 6 represents a typical unsteady flow from an animal heart. The 
Reynolds number based on average inflow velocity is 625. 

To get time-accurate solution for this problem using INS3D-UP, the solver has to be 
subiterated at each time level until a reasonable divergence free velocity is attained. In 
the present case, the number of subiterations was limited to 15. Realistically, if more 
subiteration is required, decreasing time step size while maintaining low number of 
subiteration is generally more economical. As shown in Table 1, with 200 time steps in 
each pulsating cycle, the solution still carries relatively large mass conservation error 
of up to 5 percent. Since the artificial compressibility does not require divergence free 
condition, the computational procedure produced converging solution. When the total 
number of time step per cycle is increased to 500, the resulting unsteady solution 
maintained mass conservation within a reasonable accuracy. The total time required 
for this case was 10.4 hours of Cray C-90 CPU time. 

The same case was run using INS3D-FS next. When the time step size is large with 
the total number of time step of 200, the projection error was too large. The velocity 
estimated at intermediate step and the subsequent pressure projection to the new 
velocity produced mass conservation error which failed to produce time accurate 
result. By decreasing the time step size which increases the total number of time step 
to 500, the solution maintained mass conservation and the desired time accuracy 
desired. Since the projection method does not require expensive subiterations as in 
the case of an artificial compressibility method, the computational efficiency was better 
than the result of INS3D-UP. 

The CPU time required per point per iteration varies depending on the choice of 
number of sweeps in line relaxation or iterative schemes. The comparison of 
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performance is better represented by the total time required for obtaining converged 
solutions. 


Table 1. Performance comparison for unsteady flow in curved duct flow using 
51x41x41 grid. 


Code 

Memory 

Words 

/point 

CPU/itera 

tion/point, 

nsec 

No. of 
Sub- 
iteration 



Remarks 

INS3D-UP 

40 

60 

15 

200 

4.2 

5% error in mass 
conservation 

INS3D-UP 

40 i 

60 

15 

500 

10.4 


INS3D-FS 

60 

120 

n/a 

200 

n/a 


INS3D-FS 

60 

120 

n/a 

500 

1.4 



In Figure 6, the results from these two codes are compared by showing the velocity 
vectors in mid-plane of the duct. This figure does not provide a good quantitative 
comparison. However, considering that the mass flow rate is conserved, the flow field 
inside the duct compares well between the two codes. 

4. Concluding Remarks 

The computational results show performance of the two incompressible flow methods. 
For steady state solutions, an artificial compressibility method combined with 
compressible flow algorithms is shown to be very efficient. For time-accurate solutions, 
on the other hand, the subiteration procedure required of an artificial compressibility 
method tends to be expensive. The pressure projection method offers a better 
alternative for unsteady flow computations. The examples presented here are only for 
Reynolds averaged Navier-Stokes type solutions, and did not address features 
required for a real fine mesh calculations like Large Eddy Simulations. 
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Figure 1 : Features of a wing tip vortex 
flowfield, and computational domain for 
wake vortex flow calculations. 




Figure 2 : Convergenge history for wake 
vortex case. 



Figure 3 : Axial progression of velocity 
magnitude along vortex coreline (INS3D- 
FS). 



y lx distance ( @ x/c - 1.447 ) 


Figure 4 : Comparison of velocity mag- 
nitude across wake vortex (INS3D-FS). 



y/c distance (x/c=1.587) 


Figure 5 : Comparison of velocity mag- 
nitudes from computed results of INS3D- 
FS and of INS3D-UP. 
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Figure 6. Velocity vectors at centerplane of 120-degree suare duct. Unsteady inflow 
profile is imposed from animal measurements. 




